{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 46,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "(20484, 100)\n"
     ]
    }
   ],
   "source": [
    "from neuromaps import datasets, images, nulls, resampling,stats\n",
    "neurosynth = datasets.fetch_annotation(source='neurosynth')\n",
    "abagen = datasets.fetch_annotation(source='abagen')\n",
    "neurosynth, abagen = resampling.resample_images(neurosynth, abagen, 'MNI152', 'fsaverage')\n",
    "rotated = nulls.alexander_bloch(neurosynth, atlas='fsaverage', density='10k',\n",
    "                                n_perm=100, seed=1234)\n",
    "print(rotated.shape)\n"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Nulls with parcellated data\n",
    "\n",
    "null模型函数也可以处理分割的数据，并通过接受一个额外的可选关键字参数来实现：parcellation。如果提供了这个参数，null函数会假定这是一个与所提供的数据处于同一空间的图像元组（左、右半球，如常）"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "一般来说，你已经有了预分割的数据，但为了演示，我们将使用nilearn获取一个表面分化（为10k fsaverage系统）"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 47,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "['description', 'labels', 'map_left', 'map_right']\n"
     ]
    }
   ],
   "source": [
    "from nilearn.datasets import fetch_atlas_surf_destrieux\n",
    "from neuromaps import datasets, images, nulls, resampling\n",
    "destrieux = fetch_atlas_surf_destrieux()\n",
    "print(sorted(destrieux))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 48,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "<class 'numpy.ndarray'>\n",
      "10242 10242\n"
     ]
    }
   ],
   "source": [
    "# 左脑和右脑都是ndarray\n",
    "print(type(destrieux['map_left']))\n",
    "print(len(destrieux['map_left']), len(destrieux['map_right']))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 37,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "76\n"
     ]
    }
   ],
   "source": [
    "print(len(destrieux['labels']))"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "不幸的是，此处提供的 Destrieux 地图集设计为左右半球具有相同的标签值。神经图中处理分块数据的函数始终假设标签 ID 在半球之间递增。此外，map_left 和 map_right 值是简单的 numpy 数组，neuromaps 更喜欢使用 GIFTI 图像。我们可以将数组转换为 GIFTI 标签图像，然后重新标记它们，以便标签在半球之间是连续的"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "请注意，我们在调用relabel_gifti()时设置了background=['Medial_wall']。这是因为，默认情况下，内侧壁的标签是42，而我们希望它被设置为0。（神经图谱函数假定0标签是背景，在大多数计算中它被省略了）"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 57,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "(<nibabel.gifti.gifti.GiftiImage object at 0x1773fce50>, <nibabel.gifti.gifti.GiftiImage object at 0x1773e4640>)\n"
     ]
    }
   ],
   "source": [
    "labels = [label.decode() for label in destrieux['labels']]\n",
    "parc_left = images.construct_shape_gii(destrieux['map_left'], labels=labels,\n",
    "                                       intent='NIFTI_INTENT_LABEL')\n",
    "parc_right = images.construct_shape_gii(destrieux['map_right'], labels=labels,\n",
    "                                        intent='NIFTI_INTENT_LABEL')\n",
    "parcellation = images.relabel_gifti((parc_left, parc_right), background=['Medial_wall'])\n",
    "print(parcellation)         "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 64,
   "metadata": {},
   "outputs": [],
   "source": [
    "from nibabel.gifti import giftiio\n",
    "import nibabel\n",
    "nibabel.save(parc_left,\"left.gii\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 65,
   "metadata": {},
   "outputs": [],
   "source": [
    "nibabel.save(parc_right,\"right.gii\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 70,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(10242,)"
      ]
     },
     "execution_count": 70,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# images.load_data(nibabel.load(\"left.gii\")).shape\n",
    "images.load_data(parcellation[0]).shape\n"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "我们可以使用这些图像，用neuromaps.parcellate.Parcellater类的一个实例来分割我们的数据"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 40,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "(148,) (148,)\n"
     ]
    }
   ],
   "source": [
    "from neuromaps import parcellate\n",
    "destrieux = parcellate.Parcellater(parcellation, 'fsaverage').fit()\n",
    "neurosynth_parc = destrieux.transform(neurosynth, 'fsaverage')\n",
    "abagen_parc = destrieux.transform(abagen, 'fsaverage')\n",
    "print(neurosynth_parc.shape, abagen_parc.shape)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 72,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([-1.3502014e-03,  5.1559869e-04,  2.3535038e-03,  2.2456176e-03,\n",
       "       -1.1681822e-03, -2.8523915e-03,  1.4707146e-03,  1.1858679e-03,\n",
       "       -3.0711740e-03, -2.8664269e-03,  2.4696981e-04,  1.0797515e-03,\n",
       "       -2.9965201e-03, -2.6258340e-03, -1.2355103e-03,  2.7996144e-04,\n",
       "       -1.3259611e-04, -3.2476193e-04,  2.8039730e-04,  1.1730352e-03,\n",
       "       -6.2648574e-04, -2.5964313e-04, -3.0956690e-03, -2.0425867e-03,\n",
       "       -2.2414706e-03,  2.3510500e-03,  3.3410650e-03,  3.2120552e-03,\n",
       "        4.1632741e-03, -1.0691602e-03, -2.2302552e-03, -2.7986723e-03,\n",
       "        2.0216980e-03,  6.2060350e-04, -2.0926539e-03,  2.5778895e-03,\n",
       "       -9.7296183e-04, -1.4641039e-03, -1.8901137e-03, -1.0817513e-03,\n",
       "        2.2054312e-03,  4.0947236e-04, -1.7860460e-03, -4.9383385e-04,\n",
       "        3.0656324e-03, -3.2452718e-04, -2.8216566e-03, -4.3447549e-04,\n",
       "        7.7097403e-04, -1.7246705e-03, -6.5903686e-04, -1.1652646e-03,\n",
       "       -8.5511559e-04, -3.0647888e-04, -2.0415140e-03,  1.3204279e-03,\n",
       "        4.9750938e-04,  9.8767364e-04,  4.0703896e-03, -4.8669701e-04,\n",
       "       -1.7747578e-03, -2.1700652e-03, -2.1597017e-03, -2.0246899e-03,\n",
       "       -1.6721196e-03, -1.5553853e-03,  3.9832033e-03,  1.3004806e-03,\n",
       "        4.5645023e-03, -2.9244693e-03, -2.9611134e-03, -1.0661321e-03,\n",
       "       -9.0695685e-04,  2.9345530e-03, -1.2019442e-03,  6.7350589e-04,\n",
       "        1.6607804e-03,  1.5998636e-03, -1.0552513e-03, -2.1546718e-03,\n",
       "        1.6965470e-03,  1.3302658e-03, -2.0533910e-03, -3.3256172e-03,\n",
       "        1.7101986e-04,  2.1666919e-03, -1.4698051e-03, -1.6065268e-03,\n",
       "       -4.9585005e-04,  6.5875711e-04, -4.8972375e-05,  1.7222873e-04,\n",
       "        2.0772300e-04,  9.2856889e-04, -9.8291936e-04, -3.7965362e-04,\n",
       "       -2.8118200e-03, -1.7371381e-03, -4.8746448e-04,  2.6583818e-03,\n",
       "        3.4461808e-03,  2.6069123e-03,  3.3012119e-03,  1.9060864e-05,\n",
       "       -1.9647649e-03, -2.6160129e-03,  1.7026311e-03,  1.1823134e-03,\n",
       "       -1.7169325e-03,  3.0768814e-03, -7.1045972e-04, -7.4340124e-04,\n",
       "       -6.0037611e-04, -2.0404351e-04,  1.2707853e-03, -9.8416913e-06,\n",
       "       -1.6580971e-03, -9.0806244e-04,  2.3246759e-03, -1.2698234e-04,\n",
       "       -2.4768854e-03, -6.2517723e-04,  9.1494381e-04, -2.0153159e-03,\n",
       "       -3.7748367e-04, -6.2883622e-04, -6.0542405e-04, -4.8354332e-04,\n",
       "        1.5256371e-03,  2.6276954e-03,  7.8373222e-04,  6.1012350e-04,\n",
       "        3.6375511e-03, -7.1868690e-04, -2.0954341e-03,  4.4566518e-04,\n",
       "       -1.9648003e-03, -2.1439223e-03, -1.1604582e-03, -9.5823972e-04,\n",
       "        3.1161495e-03,  1.9255460e-03,  4.3831845e-03, -2.9148581e-03,\n",
       "       -2.2515578e-03, -6.3403806e-04,  8.1840903e-05,  3.0593667e-03],\n",
       "      dtype=float32)"
      ]
     },
     "execution_count": 72,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "neurosynth_parc\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 44,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(<nibabel.gifti.gifti.GiftiImage at 0x17797a070>,\n",
       " <nibabel.gifti.gifti.GiftiImage at 0x177979b80>)"
      ]
     },
     "execution_count": 44,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "images.relabel_gifti((image1,image2))"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "*现在我们已经有了我们的分割数组，我们可以生成我们的空映射。我们使用与上面相同的调用，但提供额外的分割参数*"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 56,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(148,)"
      ]
     },
     "execution_count": 56,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "neurosynth_parc.shape"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 73,
   "metadata": {},
   "outputs": [],
   "source": [
    "rotated = nulls.alexander_bloch(neurosynth_parc, atlas='fsaverage', density='10k',\n",
    "                                n_perm=100, seed=1234, parcellation=parcellation)\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 79,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([[-2.2414706e-03, -2.2414706e-03,  2.3510500e-03, ...,\n",
       "         1.3004806e-03, -2.6258340e-03,  3.9832033e-03],\n",
       "       [-2.1700652e-03, -3.0956690e-03, -4.9383385e-04, ...,\n",
       "         9.8767364e-04, -6.5903686e-04, -2.0246899e-03],\n",
       "       [-3.0711740e-03, -1.1652646e-03, -2.9965201e-03, ...,\n",
       "        -1.7860460e-03,  2.3535038e-03, -2.9611134e-03],\n",
       "       ...,\n",
       "       [-1.6065268e-03, -2.0954341e-03,  1.7101986e-04, ...,\n",
       "         2.6276954e-03, -6.3403806e-04, -2.1439223e-03],\n",
       "       [-2.4768854e-03, -3.3256172e-03, -1.1604582e-03, ...,\n",
       "         6.1012350e-04,  3.6375511e-03, -1.7169325e-03],\n",
       "       [-6.2517723e-04, -9.0806244e-04,  1.9060864e-05, ...,\n",
       "         1.9060864e-05,  3.0768814e-03, -6.2517723e-04]], dtype=float32)"
      ]
     },
     "execution_count": 79,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "rotated"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 55,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(20484, 100)"
      ]
     },
     "execution_count": 55,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "rotated.shape"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "我们可以像以前一样将生成的数组传递给compare_images()的nulls参数"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 43,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "r = 0.416, p = 0.376\n"
     ]
    }
   ],
   "source": [
    "corr, pval = stats.compare_images(neurosynth_parc, abagen_parc, nulls=rotated)\n",
    "print(f'r = {corr:.3f}, p = {pval:.3f}')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "base",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.9.12 (main, Apr  5 2022, 01:52:34) \n[Clang 12.0.0 ]"
  },
  "orig_nbformat": 4,
  "vscode": {
   "interpreter": {
    "hash": "e8592850fc99dbb08317bbb8a8d3ad7abde2b4025207049431b7bafbb0116a2f"
   }
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}
